1. Field of the Invention
This invention relates to apparatuses, methods and systems for use in reservoir simulation. In particular, the invention provides methods, apparatuses and systems for more effectively and efficiently simulating fluid flow in reservoirs using an overlapping multiplicative Schwarz method of preconditioning adaptive implicit linear systems.
2. Background
Reservoir simulation often requires the numerical solution of the equations that describe the physics governing the complex behaviors of multi-component, multiphase fluid flow in natural porous media in the reservoir and other types of fluid flow elsewhere in the production system. The governing equations typically used to describe the fluid flow are based on the assumption of thermodynamic equilibrium and the principles of conservation of mass, momentum and energy, as described in Aziz, K. and Settari, A., Petroleum Reservoir Simulation, Elsevier Applied Science Publishers, London, 1979. The complexity of the physics that govern reservoir fluid flow leads to systems of coupled nonlinear partial differential equations that are not amenable to conventional analytical methods. As a result, numerical solution techniques are necessary.
A variety of mathematical models, formulations, discretization methods, and solution strategies have been developed and are associated with a grid imposed upon an area of interest in a reservoir. Detailed discussions of the problems of reservoir simulation and the equations dealing with such problems can be found, for example, in a PCT published patent application to ExxonMobil, International Publication No. WO 01/40937, incorporated herein by reference and in U.S. Pat. No. 6,662,146 B1 (the '146 patent”), incorporated herein by reference. Reservoir simulation can be used to predict production rates from reservoirs and can be used to determine appropriate improvements, such as facility changes or drilling additional wells, that can be implemented to improve production.
A grid imposed upon an area of interest in a model of a reservoir may be structured or unstructured. Such grids are comprised of cells, each cell having one or more unknown properties, but with all the cells in the grid having one common unknown variable, generally pressure. Other unknown properties may include, but are not limited to, fluid properties such as water saturation or temperature for example, or to “rock properties,” such as permeability or porosity to name a few. A cell treated as if it has only a single unknown variable (typically pressure) is called herein a “single variable cell,” or an “IMPES cell” while a cell with more than one unknown is called herein a “multi-variable cell” or an “implicit cell.”
The most popular approaches for solving the discrete form of the nonlinear equations are the FIM (fully implicit method) and IMPES (Implicit Pressure, Explicit Saturations) systems, as described by Peaceman, D., Fundamentals of Reservoir Simulation, published by Elsevier London, 1977, and Aziz, K. and Settari, A.: Petroleum Reservoir Simulation, Elsevier Applied Science Publishers, London, 1979. There are a wide variety of specific FIM and IMPES formulations, as described by Coats, K. H.: “A Note on IMPES and Some IMPES-Based Simulation Models”, SPEJ (5) No. 3, (September 2000), at 245, incorporated herein by reference.
The fully implicit method (FIM) assumes that all the variables and the coefficients that depend on these variables are treated implicitly. In a FIM system, all cells have a fixed number of unknowns, greater than one unknown, represented herein by the letter “m”. As a result, the FIM is unconditionally stable, so that one can theoretically take any time step size. At each time step, a coupled system of nonlinear algebraic equations, where there are multiple degrees of freedom (implicit variables) per cell, must be solved. The most common method to solve these nonlinear systems of equations is the Newton-Raphson scheme, which is an iterative method where the approximate solution to the nonlinear system is obtained by an iterative process of linearization, linear system solution and updating. The Newton-Raphson method provides that for f(x)=0, a solution can be found using an iteration:xn+1=xn−(f(xn)/f′ (xn))   (Eq. 1)
And, given an appropriate starting point, for an equation x=g(x), the iteration xn+1=g(xn) will converge to a root “a”, if |g′(a)|\1t1.
FIM simulations are computationally demanding. A linear system of equations with multiple implicit variables per cell arise at each Newton-Raphson iteration. The efficiency of a reservoir simulator depends, to a large extent, on the ability to solve these linear systems of equations in a robust and computationally efficient manner.
In an IMPES method, only one variable, typically pressure, is treated implicitly. All other variables, including but not limited to saturations and compositions, are treated explicitly. Moreover, the flow terms (transmissibilities) and the capillary pressures are also treated explicitly. For each cell, the conservation equations are combined to yield a pressure equation. These equations form a linear system of coupled equations, which can be solved for the implicit variable, typically pressure. After the pressure is obtained, the saturations and capillary pressures are updated explicitly. Explicit treatment of saturation (and also of transmissibility and capillary pressure) leads to conditional stability. That is, the maximum allowable time step depends strongly on the characteristics of the problem, such as the maximum allowable throughput, and/or saturation change, for any cell. When the time step size is not too restrictive, the IMPES method is extremely useful. This is because the linear system of equations has one implicit variable, usually pressure, per cell. In most practical settings, however, the stability restrictions associated with the IMPES method lead to impractically small time steps.
The adaptive implicit method (AIM) was developed in order to combine the large time step size of FIM with the low computational cost of IMPES. See Thomas, G. W. and Thurnau, D. H, “Reservoir Simulation Using an Adaptive Implicit Method,” SPEJ (October, 1983), p 759 (“Thomas and Thurnau”), incorporated herein by reference. In an AIM system, the cells of the grid may have a variable number of unknowns. The AIM method is based on the observation that in most cases, for a particular time step only a small fraction of the total number of cells in the simulation model requires FIM treatment, and that the simpler IMPES treatment is adequate for the vast majority of cells. In an AIM system, the reservoir simulator adaptively and automatically selects the appropriate level of implicitness for a variable (e.g. pressure, saturation) on a cell by cell basis. (See Thomas & Thurnau.) Rigorous stability analysis can be used to balance the time step size with the target fraction of cells having the FIM treatment. See Coats, K. H. “IMPES Stability: Selection of Stable Timesteps”, SPEJ (June 2003), pp 181-187, incorporated herein by reference. The computer solution for AIM systems, however, can be difficult and inefficient because of the variable number of unknowns per cell. The '146 patent by Watts describes a recent attempt to address this problem. The '146 patent describes:                “a method for performing reservoir simulation by solving a mixed implicit-IMPES matrix (MIIM) equation. The MIIM equation arises from a Newton iteration of a variable implicit reservoir model. The variable implicit reservoir model comprises a plurality of cells including both implicit cells and IMPES cells. The MIIM equation includes a scalar IMPES equation for each of the IMPES cells and a set of implicit equations for each of the implicit cells.”        
(The 146 patent at Col. 8, lines 8-16.) The '146 patent presents “a method of solving an implicit linear equation,” having the form:Ax=C,   (Eq. 2)where A is a known matrix, C is a known vector and x is an unknown vector. (the '146 patent, col. 9, lines 62-66.)
The '146 patent discloses at least three ways of solving the mixed implicit-IMPES equation. The first is a “first linear solver method” which includes these steps:                “1. Construct a global IMPES pressure matrix equation from the mixed implicit-IMPES matrix equation . . . .        2. Compute the coefficients for the saturation equations (1.2.30) at the implicit cells.        3. Solve the global IMPES pressure matrix equation for . . . a single intermediate pressure at each cell in the reservoir, and compute pressure changes based on the intermediate pressures . . . and pressures . . . prevailing at the beginning of the iteration;        4. Update implicit equation residuals at the implicit cells based on the pressures changes [of step 3]. . . ;        5. At the implicit cells, solve for improved saturations . . . in saturation equations . . . derived using a constraint of total velocity conservation between cells.        6. Update implicit equation residuals at the implicit cells and at the fringe of IMPES cells that are in flow communication with the implicit cells based on the saturation solutions obtained in step 5.        7. Determine if a convergence condition is satisfied.”        
(The '146 Patent at col. 22, line 59 through col. 23, line 28.) See also FIG. 6A of the '146 Patent and related discussion therein. The practitioner repeats steps 2-6 until convergence. (The '146 Patent at col. 23 lines 26-28.)
The second method of the '146 patent uses the following steps:                1. Construct a global IMPES pressure matrix equation from the mixed implicit-IMPES matrix equation . . . .        2. Solve the global IMPES pressure matrix equation for intermediate pressures . . . i.e. a single intermediate pressure at each cell in the reservoir, and compute pressure changes . . . based on the intermediate pressures . . . and pressures . . . prevailing at the beginning of the iteration.        3. Compute implicit equation residuals at the implicit cells based on the pressures changes . . . computed in step 2.        4. At the implicit cells, solve for improved saturations . . . and second intermediate pressures . . . by performing one or more iterations with a selected preconditioner . . . .        5. Update implicit equation residuals at the implicit cells and at the fringe of IMPES cells that are in flow communication with the implicit cells based on the improved saturations and second intermediate pressures obtained in step 4.        6. Determine if a convergence condition is satisfied.        
(The '146 Patent at col. 24, lines 19-52.) See also FIG. 6B of the '146 Patent and related discussion therein. The practitioner repeats steps 2-6 until convergence. (The '146 Patent at col. 24, lines 51-52.)
The third method of the '146 patent “effectively requires an unstructured implicit equation solver.” Steps 1-3 and 6 of the third method of the '146 patent are the same as in the second method. Steps 4 (“4.sup.II) and 5 (“5.sup.II”) of the third method are:                4.sup.II. Solve for saturations S.subj.sup.n+⅔ and pressures P.subj.sup.n+⅔ at the implicit cells while holding fixed the pressures in the surrounding fringe (of IMPES cells) to the values . . . determined during the IMPES pressure solution. Any method can be used to generate the solutions for saturations . . . and pressures . . . , but it must be able to deal with the unstructured form of the implicit cell equations.        5.sup.II. Update residuals in the fringe of IMPES cells. Since the implicit equations have been solved, their residuals will satisfy the convergence criteria.(The '146 Patent at col. 25, lines 40-57.) See also FIG. 6C of the '146 Patent and related discussion therein.        
However, the approaches to solving this problem available today have significant disadvantages. For example, it is well known that reservoir simulation matrices can best be approximately factored if they are ordered cell by cell. This means that all operations in the factorization and preconditioning step are performed as sub-matrix operations. These sub-matrices have relatively small dimensions. For FIM matrices with m unknowns per block, these sub-matrix blocks are all of order m×m. Thus for optimum efficiency the short length software loops for performing these m×m block operations may be unrolled, in order to significantly increase the overall speed of the matrix solution. By contrast with the m×m sub-matrix blocks of a FIM system, this unrolling approach for AIM systems cannot be efficiently implemented due to the varying dimensions of the sub-matrix blocks. The '146 patent does not appear to disclose an efficient unrolling technique for FIM or AIM systems.
In addition, methods without spatial overlapping can entail a slow convergence. Such methods require significant computing resources, which increases cost and lowers efficiency. The 146 patent does not appear to disclose spatial overlapping.
The present invention includes use of a global pressure matrix solution combined with efficient unrolling of loops, including spatial overlapping, in the approximate factorization of the multivariable sub-system, as a preconditioning for AIM systems.